%% Endpoint Strain analysis for Cell and Corona Data

% Rewriting strain calculation so that each time's strain is relative to
% the first timepoint 

% This will be more informative than incremental strain and cumulative
% strain, especially for oscillatory extruding cells


function [CellStrain,CellEllipse] = CellStrainCalc_endpoint(CellEllipse, tlast,cells)
%[CoronaStrain_end,CoronaEllipse] = CellStrainCalc(CoronaData)
count=0;
tic;

if (isa(CellEllipse,'table'))
    CellEllipse = table2struct(CellEllipse);
end

p = cells; 


%%Strain Calculations:
size1=tlast*length(p);
        CellStrain(size1).time=0;
        CellStrain(size1).CellID=0;
        CellStrain(size1).Exx=0;
        CellStrain(size1).Exy=0;
        CellStrain(size1).Eyx=0;
        CellStrain(size1).Eyy=0;
        CellStrain(size1).EA=0;
        B=zeros(length(p),5);
        B0=zeros(length(p),5);
        EArea=zeros(length(p));
for jj=2:tlast
    for ii=1:length(p)
        %fprintf('jj: %d p: %d .\n',jj,p(ii));
        f=[CellEllipse([CellEllipse.SliceNumA]==jj).CellID];
        g=[CellEllipse([CellEllipse.SliceNumA]==1).CellID]; 
        f = unique(f);
        g = unique(g); 
       % g=[CellEllipse([CellEllipse.SliceNumA]==jj-1).CellID]; 
       % the only change is changing the reference to the first TS rather
       % than the previous one, for all instances of jj-1 --> 1, or 2, or
       % wherever it starts 
        if(ismember(p(ii),f) && ismember(p(ii),g))
            B(ii,1)=CellEllipse([CellEllipse.SliceNumA]==jj & [CellEllipse.CellID]==p(ii)).CentroidX;
            B(ii,2)=CellEllipse([CellEllipse.SliceNumA]==jj & [CellEllipse.CellID]==p(ii)).CentroidY;
            B(ii,4)=CellEllipse([CellEllipse.SliceNumA]==jj & [CellEllipse.CellID]==p(ii)).major/2;
            B(ii,5)=CellEllipse([CellEllipse.SliceNumA]==jj & [CellEllipse.CellID]==p(ii)).minor/2;
            B(ii,3)=CellEllipse([CellEllipse.SliceNumA]==jj & [CellEllipse.CellID]==p(ii)).Angle*pi/180;
            Area=CellEllipse([CellEllipse.SliceNumA]==jj & [CellEllipse.CellID]==p(ii)).Area;
            B0(ii,1)=CellEllipse([CellEllipse.SliceNumA]==1 & [CellEllipse.CellID]==p(ii)).CentroidX;
            B0(ii,2)=CellEllipse([CellEllipse.SliceNumA]==1 & [CellEllipse.CellID]==p(ii)).CentroidY;
            B0(ii,4)=CellEllipse([CellEllipse.SliceNumA]==1 & [CellEllipse.CellID]==p(ii)).major/2;
            B0(ii,5)=CellEllipse([CellEllipse.SliceNumA]==1 & [CellEllipse.CellID]==p(ii)).minor/2;
            B0(ii,3)=CellEllipse([CellEllipse.SliceNumA]==1 & [CellEllipse.CellID]==p(ii)).Angle*pi/180;
            Area0=CellEllipse([CellEllipse.SliceNumA]==1 & [CellEllipse.CellID]==p(ii)).Area;
            EArea(ii)=(Area-Area0)/Area0;
        else
            B(ii,1)=0;
            B(ii,2)=0;
            B(ii,4)=0;
            B(ii,5)=0;
            B(ii,3)=0;

            B0(ii,1)=0;
            B0(ii,2)=0;
            B0(ii,4)=0;
            B0(ii,5)=0;
            B0(ii,3)=0;
            
        end
    end
    
    F=ComputeDeformation(B, B0);
    E=ComputeStrain(F);
    for kk=1:length(p)
        count=count+1;
        CellStrain(count).time=jj;
        CellStrain(count).CellID=p(kk);
        CellStrain(count).Exx=E(kk,1,1);
        CellStrain(count).Exy=E(kk,1,2);
        CellStrain(count).Eyx=E(kk,2,1);
        CellStrain(count).Eyy=E(kk,2,2);
        CellStrain(count).EA=EArea(kk);
    end
end
toc
end